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Abstract. The main purpose of this paper is to analyze solutions to a 
fully nonlinear parabolic equation arising from the problem of optimal portfolio 
construction. We show how the problem of optimal stock to bond proportion 
in the management of pension fund portfolio can be formulated in terms of the 
solution to the Hamilton-Jacobi— Bellman equation. We analyze the solution 
from qualitative as well as quantitative point of view. We construct useful 
bounds of solution yielding estimates for the optimal value of the stock to bond 
proportion in the portfolio. Furthermore we construct asymptotic expansions of 
a solution in terms of a small model parameter. Finally, we perform sensitivity 
analysis of the optimal solution with respect to various model parameters and 
compare analytical results of this paper with the corresponding known results 
arising from time-discrete dynamic stochastic optimization model. 
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1. Introduction and problem formulation 

In this paper we are analyzing solutions to the Hamilton-Jacobi-Bellman equa- 
tion arising from stochastic dynamic programming for optimal decision between 
stock and bond investments during accumulation of pension savings. Such an op- 
timization problem often arises in optimal dynamic portfolio selection and asset 
allocation policy for an investor who is concerned about the performance of a port- 
folio relative to the performance of a given benchmark (see e.g. [THl (THl UHl HH 13 

Consider the function V{t, y), {t, y) e T), defined on a domain T) = [0, T) x (0, oo) 
and satisfying the following fully nonlinear Hamilton-Jacobi-Bellman parabolic 
partial differential equation: 

dV ( dV 1 d'^V\ 

(la) ^ + max [A^^Xy)^^ + -,BHO,t,y)— j = 0, it,y) e V, 

and the terminal condition at t = T, 

(lb) V{T,y)^U{y), y e {0,^), 
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where U = U{d) is a smooth strictly increasing concave bounded function and e is 
a smaU parameter, < e <C 1. Moreover, we suppose that the following additional 
requirements are met: 

(1) the admissible set A* = itt] C M for aU < t < T; 

(2) the function 9 6* A^{9,t,y) e K is (not necessarily strictly) concave 
in the variable and it is function increasing at 9 = It] 

(3) the function At B 9 i-^ B^{(), t, u) is strictly convex in the 9 variable and it 
is decreasing a.t 9 — It- 

Let us suppose for a moment that the function y i— > V{t,y) is an increasing and 

strictly concave function in the y variable. Then applying the first order necessary 

condition on the maximum of the function 

r)V 1 B'^V 
9^M0Ay)^ + -^B\0,t,y) — 

we obtain the following implicit equation for 9, the maximizer of the above function: 
(2) G{9,t,y) = - where G(0,t,y) = -- 



Since the requirements ©-([H]) guarantee the increase of the function G{9, t, y) in 
the 9 variable, there exists the inverse of G and thus the unique 9 = 9{t, y) such 
that 

Then the optimal value of 9 solving ([Ta|) with the terminal condition (jlbp is given 

by 

(3) r(t,y)=min{wt,^(t,y)}. 



The problem (jlap can be now treated as a fully nonlinear parabolic partial differ- 
ential equation of the form: 

(4a) ^+^(^'^'^'^'a^) = 

where 

(4b) :F{t, y, V, ^, — ) = M9*{t, y), t, y)— + -B^{9*{t, y), t, y) — 

where 9* is given by ® and 9(t.V) = (- {^i^^v)) I {v^it^v))) depends 
itself on the solution V and its derivatives. The solution is subject to the terminal 
condition V{T, y) — U{y) where V — V{t, y) for y > and {) <t <T. Furthermore, 
— > 

The application of this study to financial markets, particularly to the theory of 
the optimal portfolio construction, has a strong impact on the special choice of the 
functions A^, and B used in the original formulation of the studied problem (jlap . 
Hence let us consider 

(5) A,{9Xy)=e+\[it{Q)-lit\y. B(9,t,y) ^ at(9)y, 

where {/3t;0 < i < T}, {iit(9)\^ < t < T} and {crt(6');0 < t < T} are assumed to 
be given deterministic processes for any choice of the control parameter 9 € At- 

Moreover, if V{t, y) is strictly convex and increasing in the y variable, then with 
regard to assumptions (|2])-([3|) the monotonicity of the function 9 i-^ G{9, i, y) is 
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guaranteed. Thus the unique maximizer 9 rising up from the imphcit equation (P]) 
fulfills 

^<'-)^«-(-j;^) — '^'')^^^^ 

Here (.)' denotes the derivative with respect to the variable and the inverse 
is taken with respect to the 9 variable. 

The paper is organized as follows. In the next section we present a stochastic 
dynamic optimization model that can be used for solving the dynamic problem of 
optimal stock to bond proportion when managing the saver's pension fund portfolio. 
We recall key steps of derivation of a time discrete version of the model. Then we 
propose a time continuous version of the stochastic dynamic optimization problem. 
We take as our setting the standard continuous-time framework pioneered by Mer- 
ton, Samuelson, Bodie, Browne and others (c.f. [H [H [201 ED 13 H H ISS IH] ) ■ It 
will be shown that the intermediate value (or utility) function satisfies certain fully 
nonlinear parabolic equation of Hamilton-Jacobi-Bellman type. In Section 3, fol- 
lowing ideas of the recent paper by Abe and Ishimura \V\ , we introduce a Riccati like 
transformation in order to transform the fully nonlinear Hamilton-Jacobi-Bellman 
equation into the quasilinear parabolic equation (see also [7]). Then we derive some 
useful bounds of a solution by using a parabolic comparison principle. These par- 
abolic estimates are also used in order to provide bounds for the original variables 
- both the intermediate utility function and for the optimal stock to bond propor- 
tion. A solution to the transformed quasilinear parabolic equation is constructed 
by means of weakly nonlinear analysis, i.e. we seek a solution in the form of infinite 
power series with respect to a small parameter representing yearly percentage of 
a salary transfered to saver's pension account. We show that first three terms in 
the expansion can be explicitly found. We also provide a recurrent formula for 
calculating higher order terms in the expansion. The rest of Section 3 is devoted to 
the sensitivity analysis of the optimal stock to bond proportion with respect to var- 
ious model parameters, in particular, to the yearly percentage transfer and to the 
saver's risk aversion. In Section 4 we demonstrate results of numerical simulations 
with model parameters corresponding to the second pillar of the pension system 
in Slovakia. We compare the results obtained with help of explicit approximations 
of an optimal solution to the time continuous Hamilton-Jacobi-Bellman equation 
with those obtained by the time discrete model from Melicherci'k et al. [9l [15] 

2. Motivation and model derivation 

In this section we first recall a discrete dynamic stochastic optimization problem 
arising in optimal portfolio selection. The discrete version of this model has been 
derived by Kilianova, Melichcrcfk and the author in [9]. It was applied for solving 
a problem of construction of an optimal stock to bond proportion in pension fund 
selection for the second pillar of the Slovak pension system. In what follows, we 
recall key steps in derivation of the discrete dynamic stochastic optimization pension 
savings model due to Melichercfk et al. in [9]. In the second part of this section 
we shall generalize the model from its discrete version to a continuous one. It 
will be shown that the continuous model for solving a problem of optimal stock 
to bond proportion in pension fund selection can be reformulated in terms of a 
fully nonlinear parabolic equation also referred to as the Hamilton-Jacobi-Bellman 
equation. 
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In the discrete optimal pension fund selection model due to Melichercik et al. 
[21 [TS] , a future pensioner with the expected retirement time in T years transfers 
regularly once a year an e-part of his yearly salary with the deterministic rate 
of growth Pt to the pension fund investing in financial market with the yearly 
stochastic return r^. More precisely, we denote by Bt his yearly salary at the year 
t. Then the budget constraint equation for the total accumulated sum Yt in his 
pensioner's account reads as follows: 

Yt+i = {l + rt)Yt+eBt+i, for t = 1, 2, T - 1, Fi = eBi. 

Supposing the wage growth f3t is known, we have the relation Bt+i = (1 + j3t)Bt 
between two consecutive yearly salaries. At the time of retiring a future pensioner 
will aim at maintaining his living standards compared to the level of the last salary 
at the retirement time t ~T. Therefore the absolute value of the total saved sum 
Yt at the time of retirement T does not represent the quantity a future pensioner 
will be taking care about. More important information for him is expressed by a 
ratio of the cumulative saved sum Yt and the yearly salary BT^ i.e. yt = Yt/Bt at 
t = T . In terms of the quantity yt representing the number of yearly salaries already 
saved at time t, the budget-constraint equation can be reformulated as follows: 

yt+i^yt{l + rt){l+ Pty^ +e, for f = 1, 2, T - 1, yi^e. 

For the sake of simplicity, we assume that the investment strategy of the pension 
fund at time t is given by the proportion E [0, 1] of stocks and 1 — 6* of bonds 
and that the fund return rt is normally distributed with the mean value and 
dispersion (t^(6') for any choice of the stock to bond proportion 9. It means that 

(7) rtie)r^ N{fiti9),aU0)), i.e. rt(0) - + (7t(0)Z 

where ^ ~ A'^(0, 1) is a normally distributed random variable having the probabilis- 
tic density function f{z) = exp(— z^/2). Both ^t and ct^ depend directly on 
the choice of parameter 9 representing stock to bond proportion in the portfolio. 
It assumed to belong to the prescribed admissible set At = [lt,ut] C [0, 1] for any 
time t S [0, T]. The admissible set At is subject to governmental regulations that 
may be imposed on the stock to bond proportion in a specific time t g [0,T]. At 
each time i, the mean value and volatility of the fund return rt can be expressed 
in terms of expected values of returns and volatilities o't^'*, cr^''^ of stocks 

and bonds as follows: 
(8) 

^,t{9) = 9^,['^ -f (1 - 0)^1"), aUe) = 9Hai'^r + (1 - 9n4'^r + 29{1 - 6)^' Qt. 

where pt S [— 1, 1] is a correlation coefficient between the returns on stocks and 
bonds at time t and the time-independent values of the parameters 
and a^") are known at time t G [0,T], they follow their relevant mutually indepen- 
dent Markov processes. 

Thus the time-evolution of the number of allocated yearly salaries can be formu- 
lated by the following recurrent equation 

(9) yt+i = G\{ytMf^t)). yi=e, 

where G\[y,rt) = £ + 2/^-^ for t = 1, 2, . . . , T - 1. 

Notice that rtiS) is the only stochastic variable appearing in the recurrent definition 
of the processes for the amount yt of yearly saved salaries. Our aim is to determine 
the optimal strategy, i.e. the optimal value of the weight 9t at each time t that 
maximizes the contributor's utility from the terminal wealth allocated on their 
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pension account, and so taking into account knowledge of the saver's utility function 
[/, the problem of discrete stochastic dynamic programming can be formulated as 

(10) maxE([/(yT)) 

subject to the constraint © where the maximum in the stochastic dynamic problem 
is taken over all non-anticipative strategies, time sequences of 9t stocks proportions, 
S = {{t,9t) I t = 1, . . . ,r}. Therefore the optimal strategy of the problem (fTO|) 
subject to ([9]) is the solution to the Bellman equation 

jU{y), t = T, 

(11) W{t,y)^^ Ez{W{t+l,Fl{e,y,Z))), t ^ T - 1, ...,2,1, 

\ SG At 

where F,\0,y,z) = G]iy, ^,ti9) + at{9)z). 

In this paper the major object of our study is the continuous version of the 
discrete model proposed above. Instead of time intervals [t,t + l],t = 1, ...,T — 1, 
representing discrete saving periods we shall assume that the proportion of the 
size ST of saving deposits is transfered to the saver account on short time intervals 
[0, t], [t, 2t], [T —T,T], where < r ^ 1 is a small time increment. The increase 
of the saver's account at time t + t can be therefore expressed as 

yt+r = F[{e,yt, Z), where Z - iV(0, 1), and 



(12) F;{e,yt,z)=ytexp( I ^,,{6) - f^t - -aUO)] t + at{0)z^) + 



- ST 



for < t < T. In the above expression for the function we have applied Ito's 
lemma (c.f. Kwok [13]) in order to generalize the discrete version of Fl with t = 1 
for the case when t is a sufficiently variable. 

Let us denote by Vit, y) the intermediate value (utility) function at time t e [0, T] 
corresponding to the saver's wealth of ?/ > saved yearly salaries in her account. 
Making use of the integral definition of the expected value {W{t + 1, F/(6', y, Z))) 
the continuous variant of the discrete backward value function for any choice of the 
small parameter < t <C 1 takes the subsequent form 



V{T,y)^U{y), t = T, 



V{t, y) = max / V{t + r, F[{e, y, z))f{z) dz, 



Therefore for the limit t = ^ 0+ we obtain 

(V{t + dt,yt+dt)~V{t, yt) 



<t <t + T <T. 



(13) 



maxE . 

0£At \ 



dt 



yt^y 







In general, we suppose that there exist functions Ag{6,t,y) and B{9,t,y) such 
that the random process yt,t & [0, T], is driven by the following stochastic differen- 
tial equation 



(14) 



dyt = A,{et,t,yt)dt + B{et,t,yt)dWt, 



where {VFt,0 < t < T} is the Wiener process. Then, by using Ito's lemma (c.f. 
Kwok [I3j) we obtain the expression for the differential dV = V{t + dt,yt+dt) — 



6 



Z. Macova and D. Scvcovic 



V{t, yt) in the form of a functfon of two independent variables t and y: 
(15) 

V{t + dt,yt+dt)-V{t,yt) 



[dV dV In d'^V 

— {t,yt) + A,{9t,t,yt) — it,yt) + -B^i9t,t,yt)-^{ 



dV 

B{e,t,yt) — {t,yt)dWt. 
dy 



dt 



Taking the conditional expectation of p5)l . the second part in the expression above 
can be omitted since stochastic variables B{0, yt)^(t, yt) and dWt are independent 
and E{dWt) = 0. Hence 

'V{t + dt,yt+dt)-V{t, yt) 



dt 



dV dV 1 n d'^V 

= ^ {t, y) +AA9t,t,y) — {t, y) + -B\0t,t,y) — {t,y). 

Letting dt 0^ , the maximum criterion in (jl3[) can be rewritten as follows 



n ^f V{t + dt,yt+dt)-Vit,yt) I 

= maxE \ Vt = y 

eeAt \ dt 

dV f dV 1 d'^V 1 

= ^) + {MO,t.y)^{t,y) + -^B\e,t^y) — it^y)}. 

In our modeling what remains is to determine the concrete form of the functions 
As{9,t,y) and B{9,t,y) driving the stochastic process for yt- Now, it follows 
from (|12p that, for an infinitesimal time increment < t — dt <ti I, we can apply 
Ito's lemma to obtain the expression for the differential dyt = yt+dt — yt- 

dyt - edt + yt {{^^(0) - (3t)dt + at{9)dWt) 

where dWt = Wt+dt - Wt = zVdi, Z - iV(0, 1), is the differential of the Wiener 
process. This way we have shown that the functions A{9, t, y) and B{9, t, y) driving 
the process (fT4| for the saver's wealth yt have the form: 

A,i9,t,y)^e+[fiti9)- Pt]y and B{9,t,y) = at{9)y. 
In summary, we have derived the following result: 



Theorem 1. The intermediate utility function V = V{t,y) satisfies the following 
fully nonlinear partial differential Hamilton-Jacobi-Bellman equation: 

and the terminal condition V{T,y) — U{y) fory > where A^{9,t,y) — e+[iit{9) — 
(3t]y and B{9,t,y) ^at{9)y. 



Henceforth, we shall make the following structural assumption on bond and stock 
average yields and their standard deviations: 



(H) 



:^ ) - gta['^] > 0, a, := [a['^]^ + [af)]^ - 2gta['^a['^ > bt , 

:^^i^)-^f)>0. 

The last condition at > bt is equivalent to the inequality crp'' — Qt^yf^ > whereas 
the first condition for bt can be reformulated as the inequality a^''"' — Qt(^[^^ > 0. 
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Notice that these assumptions are genericahy fulfilled in stable financial markets 
since expected stock returns should outperform bond returns and the correlation 
Qt between stock and bond returns is negative in typical market situations (c.f. 
[16l[9l[T5]) and the discussion on model parameters in Section 4. 



2.1. Hamilton Jacobi Bellman equation for a super-optimal solution. 

Unfortunately, due to the restriction 6 E [0, 1], the Hamilton-Jacobi-Bellman equa- 
tion (|16p illustrates a difficult problem that cannot be, in general, solved analyt- 
ically. Nevertheless, as we are approaching the problem of optimal pension fund 
portfolio construction we may relax the bound < I and allow the control param- 
eter 6 to vary over all nonnegative numbers, i.e. 9 G At = [0,oo). Taking 9 > 1 
would correspond to the situation when we allow for the so-called short positions 
in bonds because 1 — < in that case. 

Instead of the Hamilton-Jacobi-Bellman equation (|16p we shall consider a mod- 
ified problem in which we allow the control parameter 9 to vary over positive real 
numbers, i.e. 

(17) ^{t,y)+ {M9,t,y) — {t,y) + -B^{e,t,y) — {t,y)}=0. 

at 0e[o,oo) L ay 2 oy^ J 

Under the assumption ^(i, y) > 0, ^^(i, y) < 0, the hypothesis (H) and taking 

into account the definitions ^ of fit{9) and at{9), the unique solution 9{t,y) to the 
implicit relationship ^ is given by 



(18a) 9{t,y)^^- 



bt A fit %it^y) 



Hence the Hamilton-Jacobi-Bellman equation p?]) can be rewritten as follows: 

=^{t, y)+[e + yiiif^ - A + ^ A^,)] ^{t, y) 

^^^^ ^ 1 r Wi2r (s)i2n 2, 2 , 1 (W Mi^l^ 

2at dy^ 2 at %^{t,y) 

In what follows, we shall compare solutions V{t, y) of the original Hamilton- 
Jacobi-Bellman equation (fT6|l and the modified equation (fT9| . 



Theorem 2. Let At, At C M fee two admissible sets such that At C At for any 
time t G [0,T]. Let V{t,y) and V{t,y) be solutions to the corresponding Hamilton- 
Jacobi-Bellman equations with admissible sets At and At, i.e. 

(20) -g^it, y) + max [mO, t, y) — {t, y) + -B\9, y) — (t, y)] = 0, 

dV f dV 1 d'^V ^ 

(21) 7^(i,y) + ^&^\A,{9,t,y) — [t,y) + -B''{9,t,y) — (t,y) = 0, 
ot eeAt I- oy 2 dy^ ) 

fort e [0, r), y > 0, and satisfying the same terminal condition V{T, y) — V{T, y) = 
U{y) for y > 0. Then the solution V of equation k21\) is super- optimal for equation 
I.e. 

dV ( dV 1 d'^V ^ 

— (t, y) + max [A{9, t, y) — {t, y) + -B\9, t, y)^{t, y)] <0. 

Moreover, V{t, y) < V(t, y) for any t £ [0, T], y > 0. 
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Proof. The proof of the first part is rather simple and takes into account the in- 
equality 

max \A,dyV + ^S^a^yl < max \ A.dyV + Ib^OIv] = -dtV, 

because At C A(. The second part easily follows from the parabolic comparison 
principle. Indeed, let 9*{t,y) G At C At be the optimal solution to (|20p . i.e. 9* is 
the argument of the maximum operator in (I20|) . Hence 

— {t,y) + {M9*{t,y),t,y) — {t,y) + -B'{0*{t,y),t,y) — {t,y)}=O, 

BV r BV 1 B'^V ^ 

— {t,y) + [A,i9*it,y),t,y) — {t,y) + -B\9*{t,y),t,y)^it,y)] < 0. 

Therefore, applying the parabolic comparison principle (see e.g. [23) 1 we conclude 
the inequality V{t, y) < V{t, y) for any t G [0, T] and y > 0, as claimed. □ 

The above theorem enables us to refer to a solution V to the modified Hamilton- 
Jacobi-Bellman equation pT]) to as a super-optimal solution to the original equation 

m- 

3. Transformation to a quasi-linear parabolic equation 

In what follows, we shall simplify our model by assuming all the model parame- 
ters to be constant with respect to time t G [0, T], i.e. fj,[''^ ~ ^^^\^[^'^ ~ fi^'^\al''^ — 
a^^^ , cTj^'' = cr^'*^ , Qt = Q, Pt = /3- Consequently, at = a,bt — b, A/it = A/i. 

Following ideas borrowed from the recent paper by Abe and Ishimura [I] we 
introduce the Riccati-like transformation 



and the auxiliary function g{t, y) = e + ^'^^^ vit y) ~ \'^^y'^v{^^ y) where 



(22) a^^i^^^-|3+-^^l and c ^ a^^^ a^"^ ^ ' ^ ^'^ 

a 

Then equation (jl9p can be transformed to the form 

dV dV 
— {t,y)+g{t,y) — {t,y)=0. 

In terms of the transformed function if, 

(23) ^{t,y)^^(^^it,y)-^it,y)git,y) 
Moreover, the unique solution p7p to the implicit relationship ^ is given by 

(24) Oit,y)^'- + ^ 

a y(p{t,y) 

Let us introduce the following change of independent variables 

a; = ln?;, s = T — t where < s < a; G K, 
and the transformation: 

(25) tlj{s,x) — ^yip{t,y) where 7 



A/i 
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Then the original HJB equation stated for the intermediate utihty function 
V{t,y) can be reformulated for the function il;{s,x) as follows: 



9^ 



d 



(2^^) = Yd-. 



1 — ' 

dx 



for s e (0,T], X e K. The solution -tp is subject to the initial condition 
(26b) V(0,a:) ^ -7 J e\ for a; e 



Therefore the optimal 9 arising from the implicit equation ([2|) and originally 
given by (|24p . now expressed in terms of new variables (s, x) and ^p{s, x) takes the 
subsequent form 



(27) 



a 



7A/i 



aip{T - t,liiy) ' 



Remark 1. It must be remarked that the HJB equation (j26ap - (j26bp is not appli- 
cable to the Life-cycle model (see e.g. [2 ). Notice that in the Life-cycle model, the 
stock to bond ratio 9 is designed in a way it depends only on the age a of a future 
pensioner. A typical choice for 9 in the Life cycle model is 9 = 1 — a/100. Therefore 
such a stock to bond ratio 9 is independent of y variable. Taking into account (j27p . 
we obtain 

(28) ^(s, x) =0, for xeR and s G [0, T], 

i.e. ip is constant in the x- variable. Then equation (j26ap can be reworded to: 

(29) ^ = se-^^p, for xeR and s e [0, T]. 

OS 

Clearly, the above equality is impossible as the function -0 depends on s only. As 
a consequence, the Life- cycle model can not be described by the dynamic stochastic 
optimization model. 

3.1. The constant relative risk aversion (CRRA) utility function. In this 
part we discuss a suitable choice of the utility function U. We must emphasize that 
the utility function may vary across investors as it represents their attitude to risk. 
According to Arrow and Pratt the attitude to risk can be expressed in terms of 
the so-called coefficient of relative risk aversion defined as C(y) — —yU"{y)/U'{y). 
A constant relative risk aversion C{y) = d > for every y > would imply that 
an investor tends to hold a constant proportion of his wealth in any class of risky 
assets as the wealth varies. The reader is refereed to a vast economic literature 
addressing the problem of a proper choice of investor's utility function (see e.g. 
Friend & Blume [6], Pratt ^ and Young [26 ). 

In the case of a constant relative risk aversion C{y) = d > an increasing utility 
function U is uniquely (up to an multiplicative and additive constant) given by 

(30) U{y) = -y^-'^ if d > 1 , U{y) = \n{y) if d = 1 , U{y) ^ y^-'^ if d < 1 . 

The coefficient d of relative risk aversion plays an important role in many fields of 
theoretical economics. There is a wide consensus that the value should be less than 
10 (see e.g Mehra and Prescott |14j). In our numerical experiments we considered 
values of d close to 9. But it could be also lower for lower equity premium. It is 
worth to note that the CRRA function is a smooth, increasing and strictly concave 
function for y > 0. 
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For the purpose of our forthcoming analysis we consider the utihty function U{y) 
(i.e. the terminal condition for ([Tal) of the form 



(31) U{y) = -y'^-'^ where d > 1. 

The function U is a smooth strictly increasing concave function. Now it should be 
obvious that the power like behavior of the utility function U{y) = —y^~'^ leads to 
the constant initial condition (j26bp . i.e. 

(32) V(0, x) = "fd, for any x eR. 

3.2. Construction of appropriate sub- and super-solutions. In this part we 
shall derive effective lower and upper bounds of a solution ip to equation (|26ap . We 
restrict our attention to the case when the initial condition ip(0,x) is constant, i.e. 
the function U is the CRRA utility function of the form U{y) = —y^^"^ for some 
d > 0. 

Equation ()26ap represents a fully nonlinear parabolic equation of the form 

(33) al^^^^'^'^'a^'^^)- 

Notice that the right hand side of ([55]) is a strictly parabolic operator such that 

dH 1 

^(i,y,V',P,«) = y(l + ^)>y>0. 

In what follows, we shall construct positive sub- and super-solutions to the fully 
nonlinear parabolic equation p3p . The idea behind the construction of suitable 
sub- and super-solution is rather simple and it takes into account the form of the 
terminal condition U{y) for the intermediate utility function V{t, y) sXt = T. With 
regard to the presence of an advective term in the equation for the function V it 
is therefore reasonable to compare V{t, y) and the translated terminal function 
U{y + <;{T — t)) where is a positive function to be determined later. In terms of 
the transformed function ip a suitable candidate for a sub- or super-solution to 
can be therefore sought in the form: 

where ^(s) > 0, ■^(0) = 0, is a smooth function to be determined later. Assuming 
U{y) = —y^^'^ we have 

'ip{s,x) = -fd- 



l + ^(s)e-^' 
Next we calculate 

dip <;'{s)e^^ dt/j <;{s)e~^ 

" (l + ^(s)e-^)2' (i + ^(s)e-^)2' 

where (.)' stands for the derivative with respect to s. Now it is an easy calculus to 
verify the following identity: 

d_ 

dx 



1 + — 

dx 



■0// (1 + ^(s)e-^)3' 

Denote by ^ = 2a/c^,B = = /S.ii/{cy/a) — I/7. Then the lower order term 



in the right hand side Ti of equation (|35|) has the form 

^ {[A + sBe-^ - C^] i>) jdj—^^^-—^ f A + eB~ 2C^d- '^'^ 



dx^^ ^ ' (1 -h^(s)e-=^)2 \ l-h<?(s)e 
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Hence, 

(1 + <;(s)e-^)2 2?^e 



"fde~^ 1 + <,{s)e'^ 1 + i;(s)e~^ 

Clearly, the following two simple inequalities hold: 

0< — —- 0< <? 

~l + <r(s)e-^~ ~l + <;(s)e-^~ 

for any <; > and a; e M. Let £(s) and ?(s), s > 0, be solutions to the linear ODEs: 
= T ((-^ - 2C7d)i(s) - eB) , £(0) = 0, 

(35) 

-^'{s)^^{{2 + A)^{s)-eB), ^(0) = 0. 

Then it is a straightforward calculus to verify that both £(s) as well as ^"(5) are 
nonnegative and the functions 

1p(s,x) = "T^^^r , lb{s,x) — "T^^r , 

are sub- and super-solutions to the strictly parabolic nonlinear equation ([33]), i.e. 

ds± < H{t,y,i^, -^), dsij > H{t,y,^p, —, — ), 

satisfying the same constant initial condition 'ip{0, x) — ip{0, x) — for any a; e M. 

Applying the parabolic comparison principle for strongly parabolic equations 
(see e.g. |23) ) we deduce the following comparison result: 

Theorem 3. The solution ipl-s^x) to the fully nonlinear parabolic equation i33\) 
satisfies the following inequalities: 

< < ^(s,x) < < 00 

for any s G (0,T) and x G M where the functions £(s), ?(s), s > 0, are the unique 
solutions to the ODEs i35\) . i.e. 

(36) ^(s) = e(l - exp(-As))/A, £(s) = £(1 - cxp(-As))/A, 

where X = a — c^d and X~ a + , resp. 

Now, taking into account the relationships ([M)) . (P5| and integrating the above 
inequalities for the function ■0(s,a;) with respect to x we obtain the bound for the 
intermediate utility function V{t,y). 

Theorem 4. The intermediate utility function V{t, y) satisfies the following in- 
equalities: 

-{y + ^T - t))!-'^ < V{t, y)<-{y + ,{T - t)^-'' 

for any t G [0,T] and y > where the functions £(s), ^(s), s > 0, are given by 
i36\) . Moreover, the function y i— > V{t, y) is strictly increasing and strictly concave 
function, 

dV d'^V 

— {t,y)>0 and—{t,y)<Q for any t e [0,T], y > 0. 



Applying the previous theorem and using the expression ([24]) for the optimal 
value 9 we are now in a position to state useful bounds for the optimal stock to 
bond proportion in the optimal portfolio. 
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Theorem 5. The optimal value 6{t, y) describing the optimal stock to bond pro- 
portion in the optimal portfolio satisfies the inequalities: 

,37) ^ + ^fl + 3I^U»,,„)a + ^fl.i(^-') 



ad \ y J a ad \ y 

Moreover, as a consequence of the hypothesis (H), we have 0{t,y) > for any 
y > and t € [0, T] . If in addition, the coefficient of the relative risk aversion 
satisfies d > Afj,/(a — b) then for the terminal value 9{T,y) we have 9{T,y) = 
b/a + Afj,/{ad) < 1. 

It is worth to note that the condition d > A^/{a — b) can be expressed in terms 
of average stock/bond returns and their volatihties as follows: 

ni^) _ nib) 

(38) d> ^ ^ 



Cr('*)(cr(^) _ gcrib))- 

In the concrete example of the Slovak fully funded pension fund system (see Section 
4) the above constraint reads as c? > 1.78. It means that it is fulfilled in typical 
market data situations for an individual saver having the coefficient of relative risk 
aversion d greater than 1.78. 

In order to construct a solution ■0 to the problem (|26ap let us rewrite ■ip{s,x) in 
terms of the asymptotic series with respect to the small parameter e as follows: 

oo 

(39) V(s,^) = ^e"V^n(s,a;). 

The parameter < £ ^ 1 can be considered as a small parameter. In practical 
applications of the dynamic stochastic accumulation model, the value of e is close to 
0.09 (Slovak pension saving system discussed in Section 4) or £ « 0.14 (Bulgarian 
pension saving system |llp. As it could be obvious from the discussion in the 
next paragraph, the parameter £ is a natural candidate for a expansion parameter 
because we know the explicit solution ifiQ^s^x) for equation (|26ap for the vanishing 
parameter £ = 0. 

3.3. No contributions the zeroth order approximation. First of all we pay 
special attention to the first term in the Taylor expansion above, ifi{){s,x). Recall 
that due to the power like character of the utility function U{y) = — t/^"'', we have 

(40) ■0o(s, x) = jd for any s G [0, T], xeR. 

Indeed, any constant function in the x variable is a solution to (|26ap with £ = 0. 
Moreover, for the power like function U, the initial condition (j26bp is also constant 
and it is equal to jd. Therefore ipo{s, x) = jd for any s G [0, T] and a; € K. 

Let us consider the limiting case where there are no defined contributions, i.e. 
£ = 0. Then the solution ijj{s, x) coincides with t/jq and this is why the solution 
is constant in time s and spatial variable x. As a consequence we obtain, for £ = 
that the optimal stock to bonds proportion is also constant, i.e. 

_ a(^)(aW^ga(-)) + i^ 

for any i e [0,r] and y > 0. This observation is in agreement with Merton and 
Samuelson's result (c.f. [T71 [2S1 [THl [T^ ) stating that the stock to bond proportion 
is constant and it depends on saver's risk aversion only. Notice that 9{t, y) g [0, 1] 
provided that the condition ([38]) is fulfilled. 
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3.4. The first order approximation. In order to roughly approximate the func- 
tion ip(s,x) for smaU enough values of the parameter e, we use both the constant 
and the linear terms corresponding to the asymptotic expansion (j39[) to get 

d-y + ei)i{s,x) + 0{e'^) as e ^ 0+ 



(41) 



where "tpi^ajx) is an unknown function to be specified. Replacing the original 
function 'ip{s,x) by its linear approximation ()41|) above in the problem ()26ap . the 
Hamilton-Jacobi-Bellman equation for the function tpii-s, x) takes the ensuing form 



(42a) ^(. 
OS 



c 



^ 1 2(5 



for any s G (0, T) and a; G M. The solution is subject to the initial condition 
(42b) Vi(0, a;) = for any x G M, 

where d'y is replaced by the constant -ipQ and, for abbreviation, 6 stands for the 
following expression: 

(42c) S = a-dc^. 

The unique solution of the Cauchy problem (|42ap can be found in a separable form: 
(43) 

e^'^'* - 1 

V'i(s, x) = $i(s)e ^, where $i(s) ^7 for any s e (0, T), a; G R. 



Thus if the higher order terms in (|4ip are omitted, the explicit approximate 
solution of the problem (|42|) can be written as 

(44) . , - ^ 



^p{s,x) = Vo + eilJi{s,x) + O(e^) = d7|l + ■ 



-"}+0(£2) 



as £ — > 0+. Expanding the optimal value 9{t,y) and using the formula p4|) we 
obtain the first order approximation of 9(t, y) in the form 



(45) 



a ad 



1 + 

y 



Since S = a — dc^ is the same constant as A entering the expression for the lower 
bound of 6 (see (|37])) we may conclude that the first order approximation of the 
optimal value of coincides with its lower bound given by (|37p and 



3.5. The second term approximation. For the reason of better approximation 
of the function ^p{s,x) for small enough values of the parameter e, now we make 
use the Taylor expansion pop up to the second order term 

(46) '(Pis,x) = -00 + £-01(5,2:) +£^^2(5,2:) +0(£^) 

as e ^ O"*". Recall that we already have computed the first two terms ipo, ipi in the 
expansion. Namely, Vo = d'y,^i{s,x) = <I>i(s)e~^, where $i(s) = ^ (e~*^ — l). 
The function ip2{s,x) is an unknown second order expansion of the function ip to 
be determined. 

Inserting the quadratic approximation (j46p of the function ip{s, x) into equation 
(|26ap and calculating all the terms of the order 0{e^) we conclude that the function 
'4'2{s, x) is a solution to the following linear parabolic equation: 

(47) 



2 r 



1 



26- 
72 



diJ2 
dx 



{s,x) 



""Us) 
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satisfying the initial condition ■(/'2(0, x) — for a; e R, where 

(48) ^^(^)-^'(^-(d^)*'(^)-2$,(.). 

The expHcit solution of the problem (|47p can be written in a closed form: 

(49) Ms,x)^e-^^ f ^^{z)exp(c\s- z)[l + ^-^]]dz. 



The integral appearing in (j49| can be explicitly computed and it can be expressed 
as a linear combination of of three exponential functions in the s variable. 

3.6. The general asymptotic series solution for n > 3. From the straightfor- 
ward analysis of the asymptotic expansion (|39p first two terms one can deduce not 
only the separability property of the terms with respect to both variables s and x 
but the exponential contribution of the variable x to the solution's nth term. Thus 
in order to determine the nth term of the asymptotic expansion of the solution to 
the problem (j26ap let us reformulate the original solution expansion (|39p in terms 
of the special asymptotic series with respect to the small parameter e 

oc oo 

(50) V(s, x) = J2 ^"Ms, x) = J2 e"*»(s)e-"^ 

Ti=0 n=0 

with the constant zero term iIjq{s,x) = $0(5) = d'y Then the general nth term of 
the solution asymptotic expansion can be determined recursively by the following 
linear non-homogeneous first-order ordinary differential equation: 

^2 



c2 / 12 \ 

<(s) - y (^n(n-l)(l + ^)--na + 2ndj$„(s) 



2 / _ 1 \ 2 "-^^ 

(51) = y ^^^^^ E '^n-k{sW^) - n^n-l{s) + E ^n-k{s)M^) 

° ^ k=l 

for n > 1 where flk (s) is represented by the recurrent formula below 

1 ""^ 1 

(52) r!„(s) = Vl]fc(s)$„_fc(s), ^0 = -^- 

$0 ^ $0 

Then by solving the recurrent differential equation (|5ip for the nth term we obtain 

/■« 1 2na \ 

(53) ^n{s)=J en(2)exp(^y[n(n-l)(l + ^)- — + 2dn](s-z)jdz 



where 



o / X n— 1 9 n— 1 

c n(n— l)'^-^ / \r-v ^ \ 



(54) ^„(z) = - ^ ^ E '^n-ki^Wz) - n$„_i(z) + -ir E *n-fc(^)^fc(^)- 

3.7. Qualitative behavior of the optimal value 9. In this section we will be 
concerned with some useful analytic properties of the value 9(t, y). We restrict our 
attention to the first order approximation of 9 given by the leading terms in ()45p . 
i.e. 



1 + -- 



(55) ^(i,y)^^ + ^_ 

a ad I y 
We will show that even this first order approximation is capable of capturing all 
interesting phenomena that are present in our dynamic stochastic optimization 
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problem for optimal choice of the stock to bond proportion in pension fund portfo- 
lios. 

It is worth to note that (1 — exp(— (5(T — t)))/5 > for any 5 e K. Hence it is 
easy to verify that 

(56) d^/^y^<^ ^"""^ di^^^y^^^ 

for any t G [0,T),y > 0. It means that the optimal stock to bonds proportion is a 
decreasing function with respect to time t as well as to the amount y > of yearly 
saved salaries. 

3.7.1. Sensitivity of the optimal value with respect to the small param- 
eter e. First we consider the dependence of the optimal value 9 on the small 
parameter e > representing the percentage of the transfer of yearly salary to 
pensioner's account. It follows from fSS]) that 

^''^ dS^'^y^-^y 5 

for any Q < t < T and y > 0. As a consequence of the above inequality we may 
deduce that, the optimal value 9 is an increasing function in e. Taking into account 
the possible application in the dynamic accumulation pension saving model, we can 
conclude that the higher percentage e of salary transferred each year to a pension 
fund would lead to higher optimal stock to bond proportion 9. 

Supposing that the percentage e represents investor's net contributing ratio, i.e. 

e = (1 — K)e 

where n are managing costs, a regular fee charged by the the pension fund man- 
agement institutions administering investor's private pension account and e stands 
for the gross salary ratio of the financial transfer. In the Slovak pension system one 
has e = 0.09 and k = 0.01, i.e. e = 0.0891. Evidently, ^ results in 

(58) J- (t, y) = -e-f < 

an ady o 

for any Q < t < T and y > 0. It means that the increase in managing costs implies 
decrease in the stock to bond proportion, as expected. 

3.7.2. Dependence of the optimal value on the saver's risk aversion. Our 

next sensitivity analysis is focused on the dependence of the optimal value 9 on the 
coefficient d measuring saver's risk aversion. Again, it follows from (|55l) that 
(59) 

where u)(s) = 1 — 1 +sdc 



dO A^ 

dd^'^y^^-^ 



1 -I- -uj(T~t) 

y 



5 \ s J s 

Since uj{0) ~ and uj'{s) = (1 — sdc^) exp(— (5s) we may conclude that uj{T — t) is 
positive, and, consequently ^(t,y) < provided that the coefficient of the saver's 
relative risk aversion satisfies 

1 _ [a(-)p + [gWp - 2ga(-).W 
^ ' -c^T- [C7i^)]-^[CT(>>)]^1 - g^)T • 

Notice that the fraction « 306 in our market data parameter settings discussed 
in the next section. In other words condition ([TO]) is fuUfilcd for typical values 
of saver's risk aversion coefficient d ~ 10. In summary, we have shown that the 
optimal stock to bond proportion 6* is a decreasing function with respect to the 
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saver's risk aversion. In other words, higher risk aversion leads to less amount of 
stocks in saver's portfolio, as expected. 

3.7.3. Sensitivity with respect to average stock returns. Another important 
sensitivity analysis is concerned with the dependence of the optimal value 9 on the 
performance of the stock part of the portfolio. More precisely, we shall study the 
dependence of 9 on the average stock return /x^*) . Such an analysis can be useful 
when the stock market is unstable and is exposed to large variations. 

Taking into account expression for 5 — a — dc^ and a — ^^^"^ — P + (see 
one can easily verify that 

1 - e"^'* / &\ h se-^' 



36 , , 1 



l + -uj{T~t) 

y 



where w(s) = 1 — ]+Afi- 

a 



^ ' ad 

Again, as — and tf'(s) = (1 — A^^s) exp(— (5s) we may conclude that uj{T — 
t) > 0, and, consequently -^j^{t,y) > provided that 1 - A/x^T > 0. The latter 
condition can be reformulated as 

, [a(-)]^ + [cTW]^-2ga(-)aW 

Also in this case the above structural condition is fulfilled because it reads as 
^(s) _ ^(b) ^ ^ iQ g^j^(j ^(s) _ ^{b) 0.0512 in the application discussed in the next 
section. Hence we may conclude that the optimal stock to bond proportion 9 is an 
increasing function with respect to the average stock return /i^''-'. 

3.7.4. Sensitivity vifith respect to the grovifth rate. Finally we discuss the 
dependence of the optimal value 9 on the growth rate /3. Again, taking into account 
expression for 5 — a — dc^ and a = pS^^^ — /? + f we obtain 

89 , , ApLS ,^ ^ , , , l-e"-^^ se-^' 

-7r;^{t,y) = — T~^[^ — tjwhere u![s) 



d(3 ^ ad y S 

Similarly as in the previous cases, as uj{0) = and uj'{s) — sexp(— (5s) > 0, we 

may conclude that lu{T — t) > 0, and, consequently ^(t,y) > 0. Therefore the 

optimal stock to bond proportion 9 is an increasing function with respect to the 
wage growth p. 

4. Slovak pension system and calibration of the model parameters 

We have tested the proposed model on the second pillar of the Slovak pension 
system. According to Slovak legislature the percentage of salary transferred each 
year to a pension fund is 9%, i.e. e — 0.09. It means that e can be considered 
as a small parameter. We have assumed the overall time period T = 40 of saving 
of an individual pensioner. The average value of the wage growth in Slovakia for 
the period of 40 years has been adopted from the paper by Kvetan et al. [T^] and 
has been estimated (in average value) as for 5% p. a., i.e. (3 = 0.05. Similarly as 
in Kilianova et al. stocks have been represented by the S&P500 Index. For 
the purpose of the comparison of results we have taken the same time period (Jan 
1996- June 2002) yielding the average stock return /i^") = 10.28% with the standard 
deviation a^^^ = 16.90%. As the modeling of bond returns is concerned we have 
considered the term structure of the zero coupon briborQ Parameters of bond 
returns /x^''-' and their volatilities tr'^^ have been taken from ^ (see also calibration 



^BRIBOR (Bratislava Interbank Oflfered Rate) is the former term structure in Slovakia till 
1.1.2009 
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Figure 1. 3D graph (left) and countour plot (right) of the opti- 
mal value 9{t, y) computed from the first order approximation of 



1 2 3 4 5 6 0123456789 




1 2 3 4 5 6 0123456789 



y y 
a) e = 0.09 b) e = 0.f4 

Figure 2. The mean wealth E(yt) (solid line) and intervals 
E(j/f) ± a{yt) (dashed lines) obtained by 10 000 Monte-Carlo sim- 
ulation of the recurrent equation ([T2|) for different values of e. 

results of BRIBOR term structures from 24: ). We considered the average yield 
Z^^'') = 5.16% on the one year bond with standard deviation cr'^''' = 0.882% p. a. The 
correlation between stock and bond returns was set to g = —0.1151. It is the same 
correlation values as in [9] . 

In Fig. [T] we present the 3D plot as well as the contour plot of the optimal stock 
to bond proportion 9{t, y) as a function of the time t € [0, T] and the level y > 
of saved yearly salaries. As for the approximate value of 9{t, y) we considered the 
first order expansion given by (|55p . We assumed the saver's risk aversion coefficient 
d = 10. We cut-off values of the function by the upper bound 1, i.e. we in fact 
plotted the cutted function {t, y) ^ min{0(i, y), 1}. 

In Fig.[2]we present the mean wealth ^{yt) (solid line) obtained by 10 000 Monte- 
Carlo simulations of random paths {yt, t = 1,...,T} calculated according to the 
recurrent equation yt+r = Ft{(^^ Vt, ■^)j where Z ~ A^(0, 1) and 6 = 6{t, yt). We take 
one year period r = 1. The function was defined as in (fT2ll . i.e. Ff{6,yt, z) = 
yt exp ((/it(0) — (it — \(tI{9)) t -\- at(0)zy/T) + er. The dashed line represent the 
mean wealth plus/minus one standard deviation of the random variable yt at time 
t e [0, T]. The simulation were obtained by using the optimal stock to bond 
proportion 9 = 9{t, yt) depending on the value of simulated yearly saved salary yt at 
time t. For the parameter e = 0.09 (Slovakian pension system) we observe that at 
the end of simulation period t = T the averaged saved salary E(j/t) « 5.2 meaning 
that the saver following the optimal strategy given by 9{t, y) has accumulated 5.2 
multiples of his last yearly salary. On the other hand, for the higher value e = 0.14 
(Bulgarian pension system) we obtain the averaged saved salary E(yT) ~ 8.1. 
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Figure 3. 3D graph (left) and countour plot (right) of the opti- 
mal value 6{t^ y) computed by the discrete method from [9J- 




y y 



Figure 4. Comparison of the function 6{t, y) computed by the 
discrete (dashed line) and continuous (solid line) model for different 
times a) t = 1 and b) t = 35. The overall error plot of their 
difference is shown in c). 

Finally, in Fig. [3] we present results that were computed for the same model 
parameters and saver's risk aversion d = 10 but using the discrete model derived 
in [9]. Notice that in the discrete model [9] we restrict the optimal stock to bond 
proportion 6 to belong to the interval At = [0,1]. We can see and graphically 
compare the results obtained by our analytic approximation solution and those of 
discrete model have the same qualitative behavior. Comparison of the function 
9{t^ y) computed by the time discrete and continuous dynamic stochastic model is 
depicted in Fig. 21 The maximal difference 0.33 is attained at < = 1 and y ~ 1.3. 

Recently, the dynamic stochastic accumulation model derived in Section 2 has 
been adopted and numerically tested for a Bulgarian pension saving system by 
Vulkov and Koleva (see iJTj). The main difference between Slovakian and Bulgarian 
system consists in different values of the parameter e. In the Slovak funded pension 
system e = 0.09 whereas e = 0.14 for the Bulgarian funded pension pillar (see 11]). 
The impact of a higher value of the parameter e on the averaged saved sum E(dy) is 
shown in Fig. [21 Furthermore, they proposed an efficient numerical approximation 
scheme for solving the fully nonlinear parabolic equation (j26ap . They showed the 
second order of convergence of their numerical scheme. The numerical method 
can handle general class of initial conditions (|26b[) leading thus to possibility of 
considering a wider family of utility functions U than CRRA class of functions. 

5. Conclusions 

We have analyzed a dynamic stochastic accumulation model for determining the 
optimal value of the stock to bond proportion in the pension saving decision. We 
showed how the problem can be formulated in terms of a solution to a fully nonlinear 
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parabolic equation. We provided useful bounds on a solution to this equation. 

Moreover, wc performed sensitivity analysis as well as numerical simulations of the 
model. By expanding a solution to the corresponding Hamilton-Jacobi-Bellman 
equation into power series we obtained a useful first order approximation that can 
be used in qualitative analysis of dependence of the optimal strategy on various 

model parameters. 
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